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Abstract 

Nonlinear control techniques by means of a software sensor that are commonly used in 
chemical engineering could be also applied to genetic regulation processes. We provide here a 
realistic formulation of this procedure by introducing an additive white Gaussian noise, which is 
usually found in experimental data. Besides, we include model errors, meaning that we assume 
we do not know the nonlinear regulation function of the process. In order to illustrate this 
procedure, we employ the Goodwin dynamics of the concentrations [B.C. Goodwin, Temporal 
Oscillations in Cells, (Academic Press, New York, 1963)] in the simple form recently applied to 
single gene systems and some operon cases [H. De Jong, J. Comp. Biol. 9, 67 (2002)], which 
involves the dynamics of the mRNA, given protein, and metabolite concentrations. Further, we 
present results for a three gene case in co-regulated sets of transcription units as they occur 
in prokaryotes. However, instead of considering their full dynamics, we use only the data 
of the metabolites and a designed software sensor. We also show, more generally, that it is 
possible to rebuild the complete set of nonmeasured concentrations despite the uncertainties 
in the regulation function or, even more, in the case of not knowing the mRNA dynamics. In 
addition, the rebuilding of concentrations is not affected by the perturbation due to the ad- 
ditive white Gaussian noise and also we managed to filter the noisy output of the biological system. 

DOI: 10.1103/PhysRevE.72.011919 PACS number(s): 87.10.+e, 05.45.-a 
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I. INTRODUCTION 



Gene expression is a complex dynamic process with intricate regulation networks all along 
its stages leading to the synthesis of proteins Currently, the most studied aspect is that 
of regulation of initiation of transcription at the DNA level. Nevertheless, the expression of 
a gene product may be regulated at several levels, from transcription to RNA elongation and 
processing, RNA translation and even as post-translational modification of protein activity. 
Control engineering is a key discipline with tremendous potential to simulate and manipulate 
the processes of gene expression. In general, the control terminology and its mathematical 
methods are poorly known to the majority of biologists. Many times the control ideas are 
simply reduced to the homeostasis concept. However, the recent launching of the lEE journal 
Systems Biology points to many promising developments from the standpoint of systems 
analysis and control theory in biological sciences. Papers like that of Yi et al 0] , in which 
the Barkai and Leibler robustness model ^ of perfect adaptation in bacterial chemotaxis is 
shown to have the property of a simple linear integral feedback control, could be considered 
as pioneering work in the field. 

We mention here two important issues. The first one is that the basic concept of state of a 
system or process could have many different empirical meanings in biology. For the particular 
case of gene expression, the meaning of a state is essentially that of a concentration. The 
typical problem in control engineering that appears to be tremendously useful in biology is 
the reconstruction of some specific regulated states under conditions of limited information. 
Moreover, equally interesting is the issue of noise filtering. It is quite well known that gene 
expression is a phenomenon with two sources of noise: one due to the inherent stochastic 
nature of the process itself and the other originating in the perturbation of the natural 
signal due to the measuring device. In the mathematical approach, the latter class of noise 
is considered as an additive contamination of the real signal and this is also our choice here. 
Both issues will form the subject of this investigation. 

Taking into account the fact that rarely one can have a sensor on every state variable, and 
some form of reconstruction from the available measured output data is needed, a software 
can be constructed using the mathematical model of the process to obtain an estimate X of 
the true state X. This estimate can then bejised as a substitute for the unknown state X. 
Ever since the original work by Luenberger p], the use of state observers has proven useful 
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in process monitoring and for many other tasks. We will call herein as observer, in the sense 
of control theory, an algorithm capable of giving a reasonable estimation of the unmeasured 
variables of a process. For this reason, it is widely used in control, estimation, and other 
engineering applications. 

Since almost all observer designs are heavily based on mathematical models, the main 
drawback is precisely the dependence of the accuracy of such models to describe the nat- 
urally occurring processes. Details such as model uncertainties and noise could affect the 
performance of the observers. Taking into account these details is always an important 
matter and should be treated carefully. Thus, we will pay special attention in this research 
to estimating unknown states of the gene expression process under the worst possible case, 
which corresponds to noisy data, modeling errors, and unknown initial conditions. These 
issues are of considerable interest and our approach is a novel contribution to this impor- 
tant biological research area. Various aspects of noisy gene regulation processes have been 
dealt with recently from both computational and experimental points of view in a number 
of interesting papers We point out that since we add the noise 6 to the output of the 
dynamic system in the form y = CX + 5 (see Eqs. F in Section IV) it seems that its origin 
is mainly extrinsic to the regulation process, even though it could be considered as a type 
of intrinsic noise with respect to the way the experiment is performed. On the other hand, 
when writing the equation in the form y = C{X + /A), where A is a vector of noisy signals, 
one can see that the observer could estimate states that are intrinsically noisy even though 
the processes are still deterministic. 



II. BRIEF ON THE BIOLOGICAL CONTEXT 



Similar to many big cities, with heavy traffic, biological cells host complicated traffic 
of biochemical signals at all levels. Like cars on a busy highway, millions of molecules 
get involved in the bulk of the cell in many life processes controlled by genes. At the 
nanometer level, clusters of molecules in the form of proteins drive the dynamics of the 
cellular network that schematically can be divided into four regulated parts: the DNA or 
genes, the transcribed RNAs, the set of interacting proteins and the metabolites. Genes 
can only affect other genes through specific proteins, as well as through some metabolic 
pathways that are regulated by proteins themselves. They act to catalyze the information 
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stored in DNA, all the way from the fundamental processes of transcription and translation 
to the final quantities of produced proteins. 

Considering the enormous complexity of multicellular organisms generated by their large 
genomes, one can nevertheless still associate at least one regulatory element to any com- 
ponent gene. Each regulatory system is then composed of two elements at the DNA level, 
the gene that encodes a transcriptional regulator, and the target in the DNA where this 
regulator binds to, and excerts its activator or repressor function in transcription. These 
loops of interactions represent a fundamental piece to understand the functioning of complex 
regulatory transcriptional and translational networks 0, 0] . For the purpose of modelling, 
it is essential to generate simple models that help to understand elementary dynamical 
components of these complex regulatory networks as molecular tools that participate in an 
important way in the machinery of cellular decisions, that is to say, in the behaviour and 
genetic program of cells. 

Many entities in cellular networks can be identified as the basic units of regulation, 
mainly distinguished by their unique roles with respect to interaction with other units. 
These basic units are: the genes, with codifying content, also described as structural genes; 
the regulatory elements that in the old literature were called regulatory genes, which are 
smaller fragments of DNA sequences (of the order of 5 to 20 nucleotides) called operator sites 
where regulatory proteins as well as the RNA polymerase bind to; the messanger RNAs or 
mRNAs which are the products of transcription and form the template for the subsequent 
production of proteins as encoded by the corresponding gene; the forms of each protein 
and protein complexes, as well as, all metabolites present in the cell, either as products 
of enzymatic reactions or internalized by transport systems. These units have associated 
values that either represent concentrations or levels of activation. These values depend on 
both the values of the units that affect them due to the aforementioned mechanisms and on 
some parameters that govern each special form of interaction. 

This gives rise to genetic regulatory systems structured by networks of regulatory interac- 
tions between DNA, RNA, proteins, and small molecules. The simplest regulatory network 
is made of only one gene that is transribed into mRNA, this mRNA is then translated into 
proteins, which can be activated or inhibited as a result of their interaction with other pro- 
teins or with specific metabolites. Transcriptional regulators are two-head structures, one 
being the domain of DNA interaction, and the other one is the so-called allosteric domain 
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that interacts with specific metabohtes. Taking together these properties of the molecular 
machinery, one can envision that a gene encodes a protein which can regulat its own activity, 
either positively or negatively, depending on its effect in enhancing or preventing the RNA 
polymerase transcriptional activity on its own gene by means of binding to an operator sites 
upstream of its own encoding gene. Upstream here meaning before the beginning of the gene 
where transcription initiates. A mathematical model of such a biological inhibitory loop has 
been discussed since a long time ago by Goodwin and recurrently occurred in the literature, 
most recently being reformulated by De Jong Although this case could look unrealistic, 
there are simple organisms, such as bacteria, where one regulatory loop may prove essential 
as recently discussed in detail by Ozbudak et al 10]. However, already at the level of two 
genes the situation gets really complicated, mostly because of the possible formation of het- 
erodimers between the repressors and other proteins around. These heterodimers are able to 
bind at the regulatory sites of the gene and therefore can affect it and lead to modifications 
of the regulatory process. 

Recent development of experimental techniques, like cDNA microarrays and oligonu- 
cleotide chips, have allowed rapid measurements of the spatiotemporal expression levels of 



genes 



HE 



13|. In addition, formal methods for the modeling and simulation of gene 



regulation processes are currently being developed in parallel to these experimental tools. 
As most genetic regulatory systems of interest involve many genes connected through inter- 
locking positive and negative feedback loops, an intuitive understanding of their dynamics 
is hard to obtain. The advantage of the formal methods is that the structure of regulatory 
systems can be described unambiguously, while predictions of their behavior can be made 
in a systematic way. 

To make the description very concrete, it is interesting to look at well-defined, i.e., quite 
simple mathematical models that we present in the next section that refers to single gene 
cases and single gene clusters (operons). The nonlinear software sensor for such cases is 
discussed in Section IV. A three-gene case is treated as an extension to regulatory gene 
networks and shows that the method of forward engineering still works for reasonably simple 
gene networks. The conclusion section comes at the end of the paper. 
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III. MATHEMATICAL MODEL FOR GENE REGULATION 



In this section, we use the very first kinetic model of a gen etic regulation process developed 



De Jong The model in its most general form is given by the following set of equations: 

X, = Ki„r(X„)-7iXi , (1) 
Xi = Ki^i_iXi_i - 7,X, , 1 < z < n . (2) 

The parameters Kin, K21, . . . , Kn,n-i are all strictly positive and represent production con- 
stants, whereas 71, . . . , 7^ are strictly positive degradation constants. These rate equations 
express a balance between the number of molecules appearing and disappearing per unit 
time. In the case of Xi, the first term is the production term involving a nonlinear nondis- 
sipative regulation function. We take this as an unknown function. On the other hand, the 
concentration X^, 1 < i < n, increases linearly with Xj_i. As well known, in order to express 
the fact that the metabolic product is a co-repressor of the gene, the regulation function 
should be a decreasing function for which most of the authors use the Hill sigmoid, the 
Heaviside and the logoid curves. The decrease of the concentrations through degradation, 
diffusion and growth dilution is taken proportional to the concentrations themselves. For 
further details of this regulation model we recommend the reader the review of De Jong ^ . 

It is to be mentioned here that bacteria have a simple mechanism for coordinating the 
regulation of genes that encode products involved in a set of related processes: these genes 
are clustered on the chromosome and are transcribed together. Most prokaryotic mRNAs are 
polycistronic (multiple genes on a single transcript) and the single promoter that initiates 
transcription of clusters is the site of regulation for expression of all genes in the cluster. The 
gene cluster and promoter, plus additional sequences that function together in regulation, 
are called operon. Operons that include two to six genes transcribed as a unit are common 
in nature (l7l |. 

The fact that two or more genes are transcribed together on one polycistronic mRNA 
implies that we have a unique mRNA production constant and consequently we also have 
one mRNA degradation constant. In addition, the polycistronic mRNA can be translated 
into one or several enzymes, resulting in the existence of just one enzyme production and 
degradation constant, respectively. The same applies for the metabolite produced through 



by Goodwin in 1963 generalized by Tyson in 1978 



Iq and most recently explained by 
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the enzyme catalysis. Thus, if the resulting metabolite has repressor activity over the 
polycistronic mRNA (as in the case of tryptophan 2), then the model given by Eqs. ()ll2p 
could also be applied to operons and therefore it has a plausible application to the study of 
prokaryotic gene regulation. 



IV. NONLINEAR SOFTWARE SENSOR 

Numerous attempts have been made to develop nonlinear observer design methods. One 
could mention the industrially popular extended Kalman filter, whose design is based on a 
local linearization of the system around a reference trajectory, restricting the validity of the 
approach to a small region in the state space Q, The first systematic approach for 
the development of a theory of nonlinear observers was proposed some time ago by Krener 
and Isidori jj^. In further research, nonlinear transformations of the coordinates have also 
been employed to put the considered nonlinear system in a suitable "observer canonical 
form", in which the observer design problem may be easily solved 0, |2^|2^. Nevertheless, 
it is well known that classical proportional observers tend to amplify the noise of on-line 
measurements, which can lead to the degradation of the observer performance. In order to 
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avoid this drawback, this observer algorithm is based on the works of Aguilar et al. 
because the proposed integral observer provides robustness against noisy measurement and 
uncertainties. We show that this new structure retains all the characteristics of the popular 
(the traditional high gain) state observers of the classical literature and furthermore provides 
additional robustness and noise filtering and thus can result in a significant improvement of 
the monitoring performances of the genetic regulation process. 

In this section, we present the design of a nonlinear software sensor in which one Xj, for 
j 6 (l,...,n), is the naturally measured state (the most easy to measure). Therefore, it 
seems logical to take Xj as the output of the system 

y = h{X) = X, . (3) 

Now, considering the constant Kin and the function r (X„) as unknown, we group them 
together in a function Q'(^). In addition, we consider that the output function h{X) is 
contaminated with a Gaussian noise. In such a case, the model given by the aforementioned 
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Eqs. (H)) and (j2I), acquires the form: 



r : 



y = CX + 5 

where is a n x 1 vector whose first entry is and all the rest are zero, i{X) is also a 

n X 1 vector of the form [—71X1, fCj — 7jXj]-^, 5 is an additive bounded measurement 

noise, and X eW\ The system is assumed to lie in a "physical subset" S C M". 

Then, the task of designing an observer for the system F is to estimate the vector of 
states X, despite of the unknown part of the nonlinear vector (which should be also 

estimated) and considering that y is measured on-line and that the system is observable. 

A particular representation of the software sensor that we describe here is provided in 

Fig.in 
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FIG. 1: Schematic representation of the software sensor, where the output of the system is the 
input of the software sensor and the outputs of the latter are the rebuilt concentrations. 



In order to provide the observer with robust properties against disturbances, Aguilar 
and collaborators j^] considered only an integral type contribution of the measured error. 
Moreover, an uncertainty estimator is introduced in the methodology of observation with 
the purpose of estimating the unknown components of the nonlinear vector Q'(X). As a 
result, the following representation of the system is proposed 



' Xq = CX + 5 
X = 5 + £(X) 

< 

= e(A) 

yo = Xo 
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that is, in the case of the model given by Eqs. (^J and (j2I) 

Xi = Ki^iXi^i - -iiXi , 1 <i <n , 

y = , 



(4) 



where Xo is the dynamical extension that allows us to integrate the noisy signal in order 
to recover a filtered signal, while Xn+i allows us to put the unknown regulation function as 
a new state. Thus, the task becomes the estimation of this new state (a standard task for 
an observer), and therefore the function Q is related to the unknown dynamics of the new 
state. At this point, X G M""*"^, and furthermore the following equation is generated 

X = AX + B + ES , 

where AX is the linear part of the previous system such that A is a matrix equivalent in 
form to a Brunovsky matrix, B = [0, . . . , 0, Q{X)]'^ and E = [1,0,..., 0]"^. 



24|. 



We will need now the following result proven in Ref. 
An asymptotic-type observer of the system H is given as follows: 

Ao = cx + ei {yo - m) 
x = ^ + £{x) + e2 (yo-yo) 

Q = 03 {yo - yo) 
yo = Xo , 

where the gain vector 9 of the observer is given by 



So 



Each entry of the matrix Sg is given by the above equation, where Sg is a n x n matrix {i 
and j run from 1 to n), and Sij are entries of a symmetric positive definite matrix that do 
not depend on i). Thus, Sij are such that Sq is a positive solution of the algebraic Riccati 
equation 



Se { A+p\ + (a+P ] Sg = C^C 



(5) 
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In all formulas, C = [1,0,. ..,0]. In the multivariable case we must create one matrix Se 
for each block corresponding to each output. It is worth mentioning that we can think 
about this observer as a 'slave' system that follows the 'master' system, which is precisely 
the real experimental system. In addition, Sg, as functional components of the gain vector, 
guarantees the accurate estimation of the observer through the convergence to zero of the 
error dynamics, i.e., the dynamics of the difference between the measured state and its 
corresponding estimated state. One can see that i!) generates an extra degree of freedom 
that can be tuned by the user such that the performance of the software sensor becomes 
satisfactory for him. 



In 



26l | it has been shown that such an observer has an exponential-type decay for any 



initial conditions. Notice that a dynamic extension is generated by considering the measured 
output of the original system as new additional dynamics with the aim to filter the noise. 
This procedure eliminates most of the noise in the new output of the system. The reason 
of the filtering effect is that the dynamic extension acts at the level of the observer as an 
integration of the output of the original system, (see the first equation of the system S and 
the error part in the equations of system H). The integration has averaging effects upon the 
noisy measured states. More exactly, the difference between the integral of the output of the 
slave part of system S and the integral of the output of the original system gives the error 
and the observer is planned in such a way that the error dynamics goes asymptotically to 
zero, which results in the recovering of both the filtered state and the unmeasured states. 



A. Particular Case 

For gene regulation processes, which are of interest to us here, we merely apply the 
aforewritten system of equations corresponding to the asymptotic observer S 

X, = K,^,r {X,) - 7iXi (6) 
X2 = ^2,1X1 - 72X2 (7) 
X3 = ^3,2X2 - 73X3 . (8) 

The pictorial representation of this system of equations is given in Fig. El 

The values of the parameters given in Table 1, without necessarily being the experimental 
values, are however consistent with the requirements of the model. 
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FIG. 2: The genetic regulatory system given by Eqs. © - (jHl) involving end-product inhibition 
according to De Jong 9]. A is an enzyme and C a repressor protein, while K and F are metabolites. 
The mathematical model, as used by De Jong and by us, takes into account experiments where 
only metabolite K is measured. 



TABLE I: Parameters of the model 



Symbol 


Meaning 


Value 






(arb. units) 




Production constant of mRNA 


0.001 


K2,l 


Production constant of protein A 


1.0 


K3,2 


Production constant of metabolite K 


1.0 


71 


Degradation constant of mRNA 


0.1 


72 


Degradation constant of protein A 


1.0 


73 


Degradation constant of metabolite K 


1.0 


■& 


Hill's threshold parameter 


1.0 



Using the structure given by the equations of S, the explicit form of the software sensor 
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time 

FIG. 3: Numerical simulation: solid lines represent the filtered states and the dotted lines represent 
the noisy measured state for the evolution in time of metabolite K concentration. Notice that the 
initial bad estimation is due to the initial conditions that have been chosen far away from the real 
ones. This behaviour could be improved with a better knowledge of the initial conditions. The 
units of the two axes are arbitrary, i.e., the model is nondimensional. 

is: 

Xi = Xi- 7iXi + 6^2 (2/0 - yo) 
X2 = K2,iXi - 72X2 + 93{yo - yo) 

X3 = i^3,2^2 - 73^3 + 6*4(2/0 - yo) 
^4 = ^5(2/0 - ^3) , 

yo = Xo . 

Notice that this dynamic structure does not involve the regulation function. 

We can solve Eq. © and for numerical purposes we choose = 2.5 and the standard 
deviation of the Gaussian noise of 0.001. Figure El shows the numerical simulation that 
illustrates the filtering effect of the software sensor over the noisy measured state. 

On the other hand, Fig. |3] shows the results of a numerical simulation, where the solid 
lines stand for the true states and the dotted lines indicate the estimates, respectively. 
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time time 
FIG. 4: Numerical simulation: solid lines represent the true states generated by the original process 
endowed with the Hill regulatory function and dotted lines represent the estimated concentrations 
provided by the software sensor without any knowledge about the regulatory function. Plot (a) rep- 
resents the evolution of mRNA concentration in time and plot (b) the variation of the concentration 
of protein A in time. The two axes have arbitrary units. 

V. THREE-GENE CIRCUIT CASE 

In this section we extend the previous results to a more complicated case that can occur 
in prokaryotic cells. We study a more elaborated system where one regulator affects different 
promoters and transcription units. The case corresponds to the coupled regulation of three 
genes in which the metabolite resulting from the translation of gene 1 becomes the substrate 
for the synthesis of the metabolite catalyzed by the enzyme translated from gene 2, and 
similarly for gene 3, but the metabolite 3 becomes the repressor of all the three genes 
involved, as shown in Fig. 

In this case the model is given by an extension of the model given by Eqs. fjH2p . That 
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inactive 
repressor 



FIG. 5: The three-gene regulatory circuit under consideration, 
results in the following system of differential equations: 



mRNAi] 


= KiR{[Met3]) - ji[mRNAi] 


-[Enz,] 


= K2[mRNAi] ~-f2[Enzi] 




= KslEnzi] - 73[Meti] - ai[Enz2] 


mRNAi] 


= KiR{Met3) - -fi[mRNA2] 


d rrn 

-[Enz,] 


= K5[mRNA2]--f5[Enz2] 


d ^ 1 

-[Met,] 


= KQ[Enz2] - 76[Met2] - a2[Enz^] 


mRNA^] 


= KjRilMet^]) - -frlniRNAs] 


d rrn 

j^[Enz,] 


= Ks[mRNAs]-^s[Enz^] 


d ^ 


= Kg[Enz3] - i9[Met3] , 



where [mRNAi], [Enzi] and [Metj] represent the concentration of mRNA, enzymes and 



metabolites for each gene respectively. We select as the measured variables the metabo- 
lites because we want to show that through the measurement of stable molecules such as 
the metabolites, it is possible to infer the concentration of unstable molecules such as the 
mRNAs. Note that the equations are coupled through the dynamics of the metabolites. 
Moreover, we will assume that the dynamics of mRNA is bounded but unknown. 
As we showed in the previous sections our new system can be written as: 



Xi 


= X2 + di 




(9) 


X2 


= A3A3 — 73A2 — 


QilAg 


1 ■\r\\ 

(10) 


X3 


= A2A4 - 72A3 




( ■\-\\ 

(11) 


X4 


= x^ 




(12) 


X5 


= 0i(X) 




(13) 








(14) 


X, 


= K^X^ — 76^7 — 


0:2X13 


(15) 


Xs 


= K^X<j — 75X3 




(16) 


X, 






(17) 


Xio 


= 02(X) 




(18) 


Xn 


= Xyi + dj, 




(19) 


X12 


= -fCgXis - 79X12 




(20) 


Xi3 


= K^XxA^ — 78^13 




(21) 


Xu 






(22) 


Xi5 


= UX) , 




(23) 



where mRNAi = X4, mRNA2 = Xg, mRNA^ = X14 , Enzi = X3, Enz2 = Xg, Enz3 = 
Xi3, Meti = X2, Met2 = Xj, Met^ = X12, di represent the noise, (pi{X) stand for the 
unknown dynamics. In adition, the previous systems can be written in the matricial forma 
as: 

X = AX + B{X) + Ed, X e 

y = CX= (CiX^ ... C„X-)^ , (24) 
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where in this case X* G M'^* is the zth partition of the state X so that X = 
[{XY, {X'^ff and J2T=i Xi = n; A = diag[A\ ...,A"'] where A' is x A, such 
that in the equation © is invertible; C = diag[Ci, . . . , Cm], where Ci = [1,0,..., 0] 
G M^- B{X) = diag[B\Xf, . . . , B^{Xff; E = diag[^i, . . . , where Ei = [l,Q,..., 0] 

According to the scheme presented in the previous section we construct an observer 
through the following system of differential equations 



X, 


— Xr, 4- f),J Xt — XA 

— ^2 t^iil^^i y^i) 






^2 


— K-, iC-, — --vo — rvi iCo -\- f)^^( - 

— x\.3v^3 [3^2 "1^8 ' 'Jl2\^l 


- X,) 


(26) 


X-, 


— T{^x. — ^^x„ A_ f)^„( X-, — XA 

— I\2^4 [2^3 1/131,^1 ^1/ 




(27) 


X4 


= A5 + yi4(^Ai — Aij 




(28) 




= ^^15(^^1 — -^1) 




(29) 




= Xj + ^^21 (^6 — Xq) 




(30) 




= KqXs - 78^7 - 02X13 + ^22 (Xe 


-Xe) 


(31) 




= 7^5X9 -75^8 + ^23(^6 -Xe) 




(32) 




= XiQ + ^^24(^6 — Xg) 




(33) 




= +^25 (Xe — Xg) 




(34) 


ill 


= X12 + 6'79Xi2 + 6*32 (Xii — Xi) 




(35) 


-^^13 


= KgXu — 73X13 + ^^33(Xii — Xii 


) 


(36) 




= Xi5 + 5134 (Xi — Xi) 




(37) 




= ^35(Xii — Xii) , 




(38) 



where 9i stand for the observer gain values. Note, that this extension is not a direct ap- 
plication of that developed by Aguilar et al. 3] in the sense that this is a extension to 
the multivariable case. In addition, the matrix Ai is equivalent to a matrix of Brunovsky 
'orm, which guarantees the existence, uniqueness and invertibility of the matrix solution Sg'' 
27 1 . (The existence and the uniqueness of Sg^ follows from the facts that — y / — A^ is of 
Hurwitz-type and that the pair (— f-/ — Ai,Cij is observable 28^). 

Figure El shows the numerical simulation of the filtering effect of the software sensor over 
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the noisy measured state in this case. On the other hand, Fig. [7| displays the results of a 
numerical simulation of the true states (solid lines) and the estimates (dotted lines). 




5 10 15 20 

time 



FIG. 6: Numerical simulation: solid lines represent the filtered states obtained from the noisy 
measured states for the evolution in time of metabolite concentrations, where a, b and c correspond 
to metabolite 1, 2, and 3, respectively. The units of the two axes are arbitrary (nondimensional 
model). 
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FIG. 7: Numerical simulation: solid lines represent the true states generated by the original process 
endowed with the Hill regulatory function and dotted lines represent the estimated concentrations 
provided by the software sensor without any knowledge about the regulatory function; a, b and 
c correspond to molecule 1, 2 and 3, respectively. Plot (a) represents the evolution of niRNAi 
concentrations in time and plot (6) the variation of the concentration of the corresponding enzymes 
in time. The axes of the graph have arbitrary units. 
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VI. CONCLUSION 



In this research, a simple software sensor was designed for a schematic gene regulation 
dynamic process involving end-product inhibition in single gene, operon and three gene 
circuit cases. This sensor effectively rebuilds the unmeasured concentrations of mRNA and 
the corresponding enzyme. Thus, the limitation of those experiments in which only the 
concentration of the catalytically synthesized metabolite is available, can be overcome by 
employing the simple software sensor applied here. This is a quite natural case if one takes 
into account that metabolites are quite stable at the molecular level. At the same time, we 
can reproduce the concentrations of the unstable molecules of mRNA. This is a difficult task 
in experiments, despite the fact that the mRNA dynamics has been partially or even totally 
unspecified. 

The same scheme philosophy to build the observer is applied to a three-gene circuit 
with the purpose to show that the software sensor concept could be in usage in a forward 
engineering approach. In this research however, we mentioned that we were able to show 
that the observer scheme designed in for the single output case works well also in a 
multiple variable case as embodied by a particular genetic circuit given in Fig. (0). The most 
stringent mathematical requirement for this extended applicability to the multiple output 
case is described below. The linear part of the dynamic system should be a matrix by blocks 
in which each of the blocks should be of Brunovsky equivalent form. In addition, each 
subsystem corresponding to a superior block depends only on the subsystem corresponding 
to the next nearest block. This is a feature similar to the property of Markoff processes. 
The Brunovsky equivalent form of the matrix blocks Ai together with the structure of the 
corresponding output vector Q generate an observable pair {Ai, Ci), giving us the capability 
to infer the internal states of the gene network through the knowledge of its external outputs. 
However, the special Brunovsky equivalent form of the blocks leads to the possible biological 
interpretation that each block of the linear part of the differential system represents only 
that contribution of the gene regulation mechanism that comes from reactions occurring in 
a cascade fashion. 

Another important issue that we tackled in this work is related to the way of adding the 
noise to the output of the dynamic system. Even though this is a typical situation from the 
standpoint of control process theory, to the best of our knowledge it has not yet been applied 
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in the biological context of gene regulation processes. We stress that this way of including 
noise effects could have both intrinsic and extrinsic interpretations and therefore assure a 
more general approach of the noise problems. For example, in phenomenological terms, 
perturbations on the cells due to the measuring devices and the experimental conditions, 
together with the noise produced by the nature of the electronic instrumentation, could be 
equally described in this way. 

In addition, this type of nonlinear observer could be used as an online filter being ro- 
bust with respect to model uncertainties, i.e., neither a known regulation function nor the 
parameter Ki^s is required. 

This work was sponsored in part by grants from the Mexican Agency Consejo Nacional 
Ciencia y Tecnologm through project 46980-R. 
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